package gsl_hypergeometric;


#_________________________________________________________________________________________


#		This is a wrapper around my hypergeometric_cdf_P_piped binary



#			uncomment out test() and run this module directly to see it in action


#			This module runs the c++ binary (once) when the module loads, 
# 				and pipes parameters in / results out every time the function
#               is called.

#			This is loads faster than running hypergeometric_cdf_P
#				each time the function is called


#			The downside is that the code 
#				(1) uses OPC::Open2 which is a bit rocket science
#				(2) depends on the OS pipes working properly
#				(3) depends on an immediate response from the c++ programme
#                        i.e. no buffering



#_________________________________________________________________________________________

require Exporter;


our @ISA = qw(Exporter);
our $VERSION = 1.00;
our @EXPORT = qw(   hypergeometric_cdf);






use strict;
use IPC::Open2;



my ($results_hyperG, $input_hyperG);


my $process_id;
BEGIN
{
	my $BINARY_PATH = "/net/cpp-group/Leo/bin";

	# handle autovivification
    $process_id = open2($results_hyperG, $input_hyperG, $BINARY_PATH."/hypergeometric_cdf_P_piped");
}


END
{
	close($input_hyperG);
	close($results_hyperG);
}

sub hypergeometric_cdf($$$$)
{
	my ($sample_hits, $sample_pop, $background_hits, $background_pop) = @_;
	die unless defined $background_pop;
	die unless defined $background_hits;
	die unless defined $sample_pop;
	die unless defined $sample_hits;
	print $input_hyperG "id\t". join ("\t", @_) . "\n";
	my $results = <$results_hyperG>;
	chomp $results;
	return (split/\t/, $results)[1,2];
}



sub test()
{
	printf("%14s %14s %17s %17s %14s %14s %14s %14s \n", 
				"sample +ves", "sample size",
				"background +ves", "background size", "sample%", "pop%", "overProb", "underProb");
	printf("%14s %14s %17s %17s %14s %14s %14s %14s \n", 
				"-----------", "-----------",
				"---------------", "---------------", "-------", "----", "--------", "---------");
	for (0..50)
	{
		my $sample_positives = rand(50);
		my $sample_size = 500;
		my $background_positives = rand(5000);
		my $background_size = 50000;
	
		printf("%14d %14d %17d %17d %14.2f %14.2f %14.5g %14.5g\n", 
			   $sample_positives, $sample_size,
				$background_positives, $background_size,
				$sample_positives *100.0 / $sample_size,
				$background_positives * 100.0 / $background_size,
				hypergeometric_cdf($sample_positives, $sample_size,
											$background_positives, $background_size))
	}
}
#	uncomment out to test
# 
# test();
1;


